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Abstract: Transport coefficients in non-conformal second-order hydrodynamics can 
be classified as either dynamical or thermodynamical. We derive Kubo formuale for 
the thermodynamical coefficients and compute them at leading perturbative order 
in a theory with general matter content. We also discuss how to approach their 
evaluation on the lattice. 
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1 Introduction 

The theory of QCD is weakly coupled at short distances or high temperatures, but 
strongly coupled at long distances or low temperatures. One of the major goals of 
both the experimental and theoretical programs in QCD has been to understand how 
quickly this transition occurs and at what energy. A major purpose of the heavy ion 
collision program was to see if weak-coupling behavior emerges at available energies. 
Similarly, lattice studies have investigated how close thermodynamic properties come 
to their weak-coupling values as a function of the temperature. 

Broadly speaking, we can divide properties of thermal QCD into two categories: 
dynamical and thermodynamical. Most of our information on dynamical properties is 
from experiment. Experiments show [1] that at available temperatures, QCD displays 
excellent fluid behavior with remarkably low viscosity [2] . This is very different from 
weak-coupling behavior [3, 4], but roughly consistent with strong-coupling behavior 
in similar theories which we can solve [5, 6]. The story for thermodynamic properties, 
where most of our information is from the lattice, is more complex. At T ~ 150 to 200 
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MeV, thermodynamic properties such as pressure, baryon susceptibility, and (4>ip) 
show strong temperature dependence and are far from their weak-coupling values 
[7]. As temperature rises most thermodynamic quantities approach weak-coupling 
behavior, but at different rates. Quark number susceptibilities come close to weak- 
coupling behavior already at a few T c [7, 8]. Cross-correlations between strange and 
light quark numbers transition change from the expected behavior in a hadron gas 
to the behavior of a weakly coupled plasma over this same temperature range. The 
pressure takes rather longer to approach weak coupling behavior [9]. 

We feel that the more dynamical and thermodynamical quantities we have avail- 
able, the more complete and nuanced a picture of the strong to weak coupling transi- 
tion we can obtain. With this in mind, we advocate investigating the so-called second 
order hydrodynamic coefficients and their coupling dependence. As we will argue, 
some of these coefficients are thermodynamical and can be computed on the lattice. 
They also have simple weak-coupling behavior. Indeed, the main goal of this paper 
will be to compute their leading-order weak-coupling behavior in a general theory. 
We will also discuss what would be involved in evaluating them on the lattice. 

In the the next section we will review second-order hydrodynamics and explain 
how some of the coefficients of this theory are thermodynamical. Since at least 
the mid-rapidity regions in high-energy heavy ion collisions deal with QCD at small 
quark-number chemical potentials, and since the lattice can only deal well with the 
case where chemical potentials vanish, we will assume vanishing quark-number chem- 
ical potentials. We also ignore magnetic fields. However, QCD is far from a conformal 
theory in the interesting energy regime, so we will not assume conformal symme- 
try. In this case there are three independent second-order hydrodynamic coefficients 
which are thermodynamical in nature [10, 12]. In the notation of Romatschke [13, 14] 
these are k, A3 and A4. We compute their values at weak coupling and vanishing 
masses in Section 3; see particularly Eq. (3.38) and Eq. (3.39). 

As we remarked, it would be interesting to evaluate these coefficients for QCD 
on the lattice. This can be done because these coefficients all have Kubo relations 
directly in terms of finite-temperature, Euclidean correlation functions of the sort 
which can be evaluated on the lattice (without the need for analytic continuation). 
This is possible precisely because these transport coefficients are thermodynamical 
in nature. Explicit expressions for the Kubo relations are found in the next section, 
see Eq. (2.10) to Eq. (2.12). In Section 4, we present a brief discussion of how these 
Kubo relations might be applied on the lattice. In particular, we discuss operator 
normalization and contact terms. Based on this discussion, we think the evaluation 
of k should be feasible with existing techniques, at least for pure-glue QCD [15]. The 
evaluation of A3 and A4 may be prohibitively difficult since computation of three-point 
functions are always much harder than two-point functions on lattice. 
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2 Hydrodynamics 



Hydrodynamics is a general theoretical framework for describing the behavior of 
fluids locally near equilibrium. It is organized as an expansion in gradients of the 
fluid properties. (For a recent review of relativistic hydrodynamics see [16, 17]). 
At lowest (zero) order in gradients, hydrodynamics is determined by equilibrium 
thermodynamics. The state of the fluid at each point is determined by the values of all 
conserved charge densities. For QCD, these are the momentum density V 1 associated 
with the stress tensor T^ v and the charge densities Q a , a = u,d, s, . . . associated with 
the conserved 4-currents J%. Assuming local equilibrium and an equation of state 
for the pressure P in terms of the energy density and charge densities, P = P(e, n a ), 
these currents can be determined in terms of the conserved charge densities. In 
practice one uses a slightly more convenient set of variables; 1 the energy density and 
flow 4- velocity, = eu M with g^u^u" = — 1 (where g^ v is the metric tensor, in flat 
space = r]^ = Diag [—1,1,1,1]) and the number densities n a = g^u^J^. In 
terms of these the stress tensor and current at lowest order are 

T""(e, u, n) = {e+P)u^u u + Pg» v = euV + PA"" , (2.1) 
J£(e,u,n) =u^n a . (2.2) 

Here A M " = g^ u + u^u u is a projection operator onto the local spatial directions. 
When T^ u and satisfy these expressions, then stress conservation, V^T^" and 
current conservation, V M J£, close and completely determine the fluid dynamics. The 
relevant quantities at this order - the pressure P and its various derivatives which 
give the entropy density s, quark number susceptibilities Xab, the speed of sound c s , 
and so forth - have been extensively studied on the lattice [7, 8]. 

At first order in gradients there are two independent terms one can add to the 
righthand side of Eq. (2.1): 

= RHS of Eq. (2.1) - r)a>* v - (A^ v V a u a , (2.3) 
^ _ £jj,a£vp (^j aU ^ + v ^ _ 2 Aq/3 v 7U 7^ . (2.4) 

Here rj, ( are the shear and bulk viscosities respectively. These new coefficients r), 
(, and the diffusion coefficient which can be added to Eq. (2.2), are all dynamical 
quantities. If the coefficients and the terms they multiply are nonzero then entropy 
increases. They can be determined, via Kubo relations, from equilibrium correlation 
functions of stress tensor or current operators, but the relations involve evaluating 
these correlation functions at nonzero frequency, which makes a direct evaluation on 
the lattice impossible and an indirect evaluation at best very challenging [18-20]. 
However there is significant progress in determining them from experiment [21]. 

1 Note that we are using the Landau-Lifshitz frame. 



- 3 - 



At second order there are a host of terms which can be added. The situation 
improves somewhat if we assume that charge densities are small, so J£ terms can 
be neglected. If in addition we assume that the theory under consideration is con- 
formally invariant, then there are 5 additional terms which must be included [14]. 
However, since we are interested in QCD at finite coupling and potentially in mak- 
ing contact with the lattice, we cannot assume conformal invariance. In this case, 
after reducing the number of terms by applying equations of motion and other inter- 
relations, there are 15 independent terms which appear at second order, which have 
been enumerated by Romatschke [13]. To write these terms explicitly it is convenient 
to introduce the vorticity tensor W, 

m „u = A „a A W3 ( Vq ^ _ V/?Mq ) ; (2.5) 

as well as the curvature tensor R^ ua/3 and Ricci tensor R^ u = R^ a ua and scalar 
R = And we will write 

RI M(ya)P = }_ R ^P ( A " A <* + A^A" - ^A™A re(J J (2.6) 

and similarly for R^ u \ That is, the indices enclosed in angle brackets are space- 
projected, symmetrized, and trace-subtracted. Using all of this notation, the possible 
second-order terms, according to Romatschke [13], are 

= Eq. (2.3) + 777V (u • Va"" + ^p*"") 

+k (R { v) - 2u a upR a{ v)^ 
+\ 1 a x ^a^ x + \ 2 a x ^ x - X 3 Q X { ^ X 
+r]T*^—^aV + A 4 V (/X In sV u) In s + 2K*u a upR a{ v)P 

(-Cr u u ■ W • u + ^a a(3 + 6(V • uf 
+e 4 V Q± In sVl In s + 6^0^ + &R + & Wi^) . (2.7) 

There are several ways to categorize these terms. Some are only relevant in 
curved *, £5 and ^6- The others are relevant in flat or curved space. (Even 

though k etc. only play a role in curved space, they mix with the other terms when we 
find Kubo relations in Eqs. (2.10-2.12), so they should generally be considered anyway 
[14].) We can also divide the terms into linear and nonlinear terms. Linear terms 
affect small fluctuations and can, for instance, influence their dispersion; nonlinear 
terms are only relevant at second order in small fluctuations about equilibrium and 
flat space. The linear terms are r n , n, k* , r n , £5, and £ 6 . The other terms, A1...4, r*, 
£1. 4 are nonlinear. 

We can also group these terms into those which are thermodynamical in nature, 
and those which are dynamical. We call a term thermodynamical if it can give a 
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nonzero contribution to T^ v when the geometry and density matrix are fully time- 
independent and the system is therefore in equilibrium. No term involving the shear 
tensor is thermodynamical because a system under shear flow is changing with 
time and is producing entropy. In nonconformal theories, the same is true of bulk flow 
V-it. However, it is completely consistent to have a system which is in equilibrium in 
a curved (but time- independent) geometry. Similarly, a time-independent but space- 
varying goo (gravitational potential) makes V^s nonzero without any departure from 
equilibrium. Similarly, it is possible (in a curved geometry) to establish persistent 
vorticity which is sustained forever. 2 The system will be fully in equilibrium in the 
presence of this vorticity. Hence, the coefficients k, k*, A 3 , A 4 , £ 3 , £ 4 , £5, £ 6 represent 
thermodynamical quantities. 

In Ref. [22] we showed how to derive Kubo relations for second-order hydro 
coefficients. There we did so only for conformal theories, but it is straightforward to 
do so for nonconformal theories as well. Doing so, we find that the Kubo relations for 
the thermodynamical coefficients can all be expressed in terms of retarded correlation 
functions evaluated directly at zero frequency. Up to powers of i, zero-frequency 
retarded correlators equal zero- frequency Euclidean correlators. In fact, we can derive 
(Kubo) relations between the thermodynamic coefficients and Euclidean correlators 
by working directly in Euclidean space. In particular, defining the Euclidean n-point 
function as 



G^ vl -^ Vn {p^ . . , Pn _ 1 , -pi-. . -p„_i) = / . . . rf 4 x„_ie 



-i(pi-Xl + ...+p n -l-Xn-l) 



X 



2 n d n In Z 



dgu 1 u 1 (x 1 )...dg IMtVn (p) 



with 



we find 



z[g, 



/IV J 



V(p ex.p{-S E [(f), g^u]} 



lim 



<)- 



case ™ 



A 



fc =0 

d 2 



lim -J^GT' xy {p,q)\ 



2 p x "qv^o dp x dq y ~ E 



Po,<?o=0 



(2.8) 
(2.9) 

(2.10) 
(2.11) 
(2.12) 



The remaining transport coefficients, including n* which appears above, are not inde- 
pendent but are determined in terms of these three via five independent conditions. 



2 For instance, consider a spacetime which is S 2 x 1Z 2 , with time in one of the flat directions. 
The fluid can spin about the equator of the S 2 and this flow will persist forever, and will therefore 
reach equilibrium. 
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Two conditions were found by Romatschke [16], by demanding that the entropy cur- 
rent have non-negative divergence. His calculation was limited to second order in 
gradients; but a treatment to third order in gradients by Bhattacharyya [10] and 
Jensen et al [11] found three more constraints on second order transport coefficients. 
For an interesting physical interpretation of these constraints, see [12]. 
The five constraints found by Bhattacharyya, in our notation 3 , are 

6 = \ (&% - * - I) . (2-14) 
9 / ^dn dn* \ An* A 4 

6 = ( 3r 5r " 2T It + 2k ' 3k ) ~ K + T + W ■ (215) 

i<?.T/dK" dn\ 3(<?.-l). . , A„ 1 / 2 <iA 3 2 A 3 



(2.16) 

(^) • < 2 ' 17 ' 

We take these constraints to determine all other coefficients in terms of k, A3, and A4. 
In Appendix A we give a detailed derivation of Eqs.(2.10-2.12), and we find Kubo 
relations for the dependent transport coefficients mentioned in Eqs. (2.13-2.17) for 
completeness. 

Euclidean correlation functions have well behaved perturbative expansions at 
finite temperature (at least at low order), therefore it should be possible to evaluate 
these correlators perturbatively in a weakly coupled theory. We present the deriva- 
tion at lowest order in a general massless theory for particles of spin zero, half and 
one in the next section. 



3 Evaluation at Weak Coupling 

To carry out the calculation of these transport coefficients, first we have to clarify 
the nature of the correlation functions that are derived by differentiating the curved 
space partition function. The definition for the n-point Green function established 
in Eq. (2.8) involves multiple derivatives acting on the energy functional. Each 

3 Labeling the coefficients of [12] with a prime, the relations between their coefficients k[, k' 2 , A3, 
X' 4 , C3; £3; ancl £4 an< l our coefficients arc: Tk[ = k, Tk' 2 — 2k — 2k*, — TA 3 = A3, c*TA 4 = A4, 
TC2 T C 3 = &, -T£,' 3 = £3 and c*T£ 4 = £ 4 . Moreover, unlike [12] our convention for RW is 

rp — 8 vp — 8 vp +v p r A — r p r A 



- 6 - 



derivative can pull down a factor of —IdL^jdg^ = T MV , giving a conventional n- 
point stress-tensor correlator; but the g^ u derivatives can also act on T a/3 factors 
pulled down by previous g a p derivatives, leading to contact terms. (In intermediate 
steps our T^ v is really the stress tensor density JgT^; the distinction is irrelevant 
in final expressions since in the end we evaluate correlators in flat space.) So in terms 
of the usual n-point stress tensor correlators, Q^>-> a P defined in Eq. (2.8) is 



G^(0,x) = (inojT^x))! + 2 



G ^,7P(o,x,y) = (T^{0)T^(x)T' p {y)) 
+2 
+4 



/ dT^(O) Ta0{ 

\ d 9i P {y) 

d 2 T^(0) 



x 



dg a f}(x)dg lp (y) 



gphv—'nijbv 



<9T^(0) 
dg a p{x) 



2 ( T^(0) 



dT a(3 { 



X) 



dg 1P (y) 



(3.1) 
(3.2) 

Qjiv—^tiv 



The terms involving derivatives of the stress tensor are called contact terms, and are 
discussed in some detail in Ref. [23]. Since dT pv (x) / dg a p(y) oc 5 4 (x — y) they have 
very simple momentum dependence. In particular, the contact term in G^' a/3 (x) is 
oc £ 4 (x); so its contribution to (k) is k- independent. Therefore it does not 

contribute to Eq. (2.10). 

Now consider the four contact terms in Eq. (3.2) and their contribution to 
Eq. (2.11). Defining 

2^—j^ = X^S\x-y), (3.3) 

we find three X-type contact terms, involving 5 4 (x), S 4 (y), and 5 4 (x — y) respectively. 
The first gives a contribution which is independent of p and so does not contribute 
to Eq. (2.11); similarly the second is independent of q and also does not contribute. 
But the third term does contribute; 



A, 



-4 lim d Pz d qz (T* t (p)T*(q)T»(-p - q)) 



-4 \\m d Pz d qz {X xt y\p + q)T x y{ 



-p 



?)) • 



(3.4) 



In order to calculate these transport coefficients for a generic field theory, we 
need to find the explicit form of both the stress tensor T pu and of the contact term 
j^/iuaf) differentiating the action 



S 



J d X (/^scalar 



"I - -^spinor C* 



vector J 



(3.5) 
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with respect to the metric. Since we only attempt a leading-order calculation here, 
it is sufficient to consider the free-theory action in curved space, 



Scalar = ^^,0^0, (3.6) 

= |e| ^ C e x c (d x + ^G ab co^ V , (3.7) 

£ vector = ^ Ffii/Fprd^^ 9 • (3-8) 

Here g^ u is the inverse of g^ u , F^ v = d fM A u — d u A^ is the field strength tensor, e™ is 
the vierbein related to g^ by i] a \ 1 e a lL e h v = g^ and \e\ = det (e^) is its determinant. 
Finally, u^ b is the spin connection and G ab = \ [7°, 7 b ] . 

Actually, more generally the scalar Lagrangian density should read 

^scalar = ^ {fdyQd^ ~ , (3.9) 

where £ is a dimensionless constant and R is the Ricci scalar introduced earlier. The 
action is conformal for the choice £ = | and is called minimally coupled if £ = [24]. 
We will consider general £, but in the end our results for k, A3 are £ independent. 

Through the reminder of this section we will compute k and A3 using Eq. (2.10) 
and Eq. (3.4), applying the action in Eq. (3.5). We also note that at the leading order 
calculation, A4 = due to conformal symmetry. Since the effect of more degrees of 
freedom Nq, N1/2 and N\ at this level is multiplicative, they can be counted in the 
final result respectively. 

3.1 Scalars 

In carrying out the variation of Eq. (3.9) with respect to we must consider the 
explicit dependence and the implicit dependence via the Ricci scalar R. The resulting 
stress tensor is 

T>» = ((1 - + ^-V^) 0« W + 2£ brg^-g^g"*) ^drf , 

(3.10) 

plus terms which vanish in flat space. 

With T^ u in hand, we can compute the scalar contribution to k. The lowest 
order diagram is shown in Figure 1. The momentum k enters at one T xy insertion 
and exits at the other, so the scalar propagators carry different momenta, p and 
p+k = q. T xy inserted between these lines obeys the Feynman rule (note directions 
of momentum flow) 

= (1 - 20 (p x q v + q x pA + 2c(p x p y + fA . (3.11) 



p^ q 
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p+k 




P 



Figure 1. Leading order scalar diagram contributing to (T xy \—k)T xy (k)) , necessary for 
evaluation of k. The crosses are T insertions, the solid lines are scalar propagators, and 
the arrows indicate the momenta flowing on lines and entering or leaving T insertions. 



Since we are only differentiating with respect to k z , we may set k x = = k y from 
the outset, in which case p x = q x and p y = q y . Therefore the £ terms cancel 4 and 
the diagram evaluates to 

K = dl(T xy (-k)T xy (k)) = -d 2 k ^ ( 2pXpV "> 2 



2 * z %P 2 {p+k) 



k=0 



-^t{^-f)Apl = ~ (3.12) 



where \ is the symmetry factor of the diagram, and the integration-summation sym- 
bol is defined as 

d 3 p 



Jp p°=2nnT 



(2vr) : 



(3.13) 



and n runs over the integers. In evaluating this and related sum-integrals we use the 
result 

4 (P 2 T (2n + l)! T 2 

% (p2)n+i 2 2 ™(n!) 2 12 ' 1 ' 

Expressions with powers of (p ) 2 in the numerator can be handled by rewriting 
ijp) 2 = p 2 — p 2 and using this relation repeatedly; for instance, ^ j^r = -j^ We 
handle p 2 x p 2 p 2 by angular averaging, (p^p^p 2 ) an gie = p 6 /105. 

Our leading-order result for k agrees with the result in Ref. [23]. The above 
result shows that the weak coupling expansion of k starts at a . 

Now we turn to the computation of A3. The first term appearing in Eq. (3.4) is 
represented by the diagram shown in Figure 2. Once again, p, q only need have non- 
vanishing z-components, but the T operators only return x, y, components, which 
simplifies the evaluation of the diagram and ensures that the result is £ independent. 

4 The same will not happen if we compute non-conformal coefficients. 
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Figure 2. Three-point correlation function {T xt (p)T yt (q)T xy (— p — q)) that contributes to 
the Kubo formula of A3; the leftmost vertex is T xy , the other vertices are T xt and T yt . 



The diagram evaluates to 



-4d p ,d q , (T xt (p)T yt (q)T xy (-p - q)) 

JL (2k x k y ){2k t k x )(2k t k y ) 

~ ~ A6pA % (k+p) 2 (k-q) 2 k 2 p 

r k 2 kl kl k 2 t 2 
= ]3 g£j^_ M = __. (3 . 15) 

For the contact term in Eq. (3.4), we need to calculate X ytxt as denned in Eq. (3.3). 
Variation of Eq. (3.10) with respect to g xt gives gives 

X ytxt = -d y <f)&*(f), (3.16) 

in flat space and for £ = 0. In addition there are terms proportional to £, but they 
again always involve the combination d a dp(j) 2 . These terms do not contribute to the 
correlation function we need for the same reason the ^-proportional terms above did 
not contribute; the incoming and outgoing momenta equal for the components which 
make up the indices of X txty . Therefore the result is again £ independent. The 
contribution from the (TX) correlator to A3 is 



- 4d Pz d qz (X°**(p + q)T*y(-p -q))= 2d Pz 8 q ^ ^ 



[2k x ky)' 



p) 2 (k+q) 2 
= - 32 ±^^ = -—. (3.17) 

This diagram is actually the same as the diagram which determines k; shifting the 
integration variable used above by p, the integral becomes the same one needed in 
evaluating k except for the overall factor of 4. Summing up Eq. (3.15) and Eq. (3.17), 
we get 

T 2 

A3 = 1 real scalar field. (3.18) 

3 12 v ; 

We follow a similar calculation for gauge and fermion fields in the next sections. 
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3.2 Gauge fields 

The gauge field stress tensor derived by g^ v variation using the gauge field action, 
Eq. (3.8), is 

T r> = F^ a p v a - I^F a/3 F Q/3 , (3.19) 

and from the above relation, we derive the Feynman rule for the vertex, 

/i, a 'wtgryvp u, b 5 ab ( (p a g^ - p^g^a) (kpgl - k"'g u/3 ) + (p u) 
p k-*~ v 

- g a p(p ■ k g^ - kpPv)) . (3.20) 

The expression for X^ a/3 is rather long, but for the case of interest, X txty , it is quite 
simple: 

x txt y = _ F xz F y z ^ 321 ^ 

leading to the Feynman rule 

X 0x0y 

fi,arw<8ryvnv,b - 5 ab [(p 

xg^iz Pzgfj,x){kyg l , z k z g V y) 4~ [x -H- y) ) . (3.22) 
p -*~ k-*~ v / 

The calculation of k and A3 then proceeds via the same diagrams as in the scalar 
case, but with these somewhat more complicated Feynman rules for the vertices, 
and with gauge propagators. Note that, because T^ u is built from field strengths, 
it applies a transverse projector onto the incoming gauge field index; contracting 
Eq. (3.20) with p^ or k v gives zero. Therefore the result is gauge parameter inde- 
pendent within covariant gauges (and all linear gauges). After significant algebra we 
find that 

T 2 

k = — — for a single color , (3.23) 
18 

while the two diagrams contributing to A3 give 

2T 2 

- 43,A {T*Xp)TvXq)T*v{- V - q)) = — (3.24) 

and 

- m Pz d qz (x^( P + q)T*y(- P -q)) = Y- ( 3 - 25 ) 

Therefore, the gauge field contribution to A3 is 

T 2 

A 3 = — for a single color. (3.26) 
3 

Since each color possesses two spin states, we need to divide these results by 2 to get 
expressions per degree of freedom. 
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3.3 Fermions 



The treatment of fermions in curved space requires the introduction of the vierbein 
(also called the frame vector or tetrad) (for a review and a treatment of their 
application to the stress tensor see Ref. [25]). The vierbein relates a local orthonormal 
coordinate system on the tangent space, with indices a and metric j]^ (which is S a j, 
in Euclidean space) to the metric, via 

Wle h v = gilv ; (3.27) 

in a sense it is the square root of the metric. The Dirac action is e^7 a V ^ip , where 
the action of V M on a spinorial object is determined by the spin connection cu^ b : 

= + -G^V (3.28) 
where G[ a b\ = \ [7a, 7ft], and the spin connection is related to the vierbein via 

< = - d v e\) - \e bv {d,el - d v e a J + l -e™ e b ° {d a e c v - d v e%)e ctl . (3.29) 

Because the action depends on the local frame components, the stress-tensor for 
fermions cannot be obtained by functional differentiation with respect to the metric 
tensor; instead one must use the more general expression T» u (x) = e u a g£^- r which 
reduces to T^ u (x) = 2 g g dZ ^ for any terms which depend only on g^ u because of 
Eq. (3.27). Applying this relation to the fermionic action, noting that = —rf h g ilv 
(since g^ is the inverse of g^ v \ alternatively, because the variation of g^g^v should 
vanish), and specializing to the non-diagonal entries in T^, after some work one 
obtains [25] 

T fW = - V v ip - 7^ + - ■ (3-30) 

The relevant Feynman rule is 

JT^Z I W + «") + ^ + <f)) , (3-31) 
leading to an expression for k, 



Tr ([(2p+k) x 7 y + {2p+k)y-i x }i>[(2p+k) x -i y + (2p+k) y 7 x ]\rf+ft}) 

k z " 



(3.32) 

where prime on the sum-integral indicates that the frequencies are (2n + 1)ttT. In 
terms of this fermionic sum-integral, the equivalent of Eq. (3.14) is 

' (P 2 Y (2n+l)! (-T 2 ) 
(p2)n+l 2 2 ™(n!) 2 24 ' K ' 
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The rest of the evaluation is straightforward, yielding 

T 2 

k = — for a single flavor (3.34) 

Since a Dirac fermion has 4 degrees of freedom, this should be divided by 4 to get 
the contribution per degree of freedom. 

Next we calculate A3 from Eq. (3.4). The three point diagram still looks like 
Figure 2, but with two terms depending on whether the fermion number follows 
or opposes the indicated momentum flow. A straightforward evaluation yields a 
contribution to A 3 of T 2 /24. 

We specialize immediately to the contact term needed in the calculation; the 
general expression is not simple. The formula for X txty is 

1 / ST ty ST ty 5T tx 5T tx \ 
X txty = 7 U'j- + ef-f- + e*£- + . (3.35) 

4 V oe ax de at de ay de at J 

Using the techniques already introduced, we find after some work that 

x txt y = Z_ (^VV + $l v V x ip - V x h v ^ ~ VVf ^) = ^-T xy . (3.36) 

The (TX) diagram contributing to A3 therefore gives 3 times the contribution we 
found for the diagram contributing to k, that is, T 2 /24. Adding these two terms, we 
find 

T 2 

A3 = — for a single flavor. (3.37) 
Once again, we need to divide by 4 to get the contribution per degree of freedom. 
3.4 Results 

Since we work, so far, at the free theory level, the result is a function only of the 
number of scalar, spinor, and vector degrees of freedom, which we will write as Nq, 
Ni , Ni. 5 Combining the results of the previous subsections, we find 



T 2 

K 



288 
T 2 



(-4/V + Ni + 8N^j + O(y/o) , (3.38) 

^ (-4iV + Ni + SiVx) + O(y/o) ■ (3-39) 

The other coefficients vanish because the theory is conformal at this order. Curiously, 
at the free level A 3 = 6n regardless of the matter content. 

We have computed only the leading, coupling-independent contributions. We 
expect the first corrections to k, A3 to arise at O(a^), and the first contributions to 
the nonconformal coefficients to be O(a). 



5 No is 1 per real scalar field; Ni is two per Weyl spinor field, or 4 per Dirac field; and is two 
per massless spin-1 field (one per spin state). For 3-flavor QCD, No = 0, Ni = 4 x 3 x 3 = 36 [4 
for a Dirac spinor, times three colors times three flavors], and Ni = 16 [2 spin states times 8 color 
combinations]. For U(N) Af=4 SYM theory, N = 6N 2 , Ni = 8N 2 , and Ni = 2N 2 . 
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4 Lattice implementation 



Here we will briefly discuss some of the challenges associated with evaluating the 
second-order coefficients on the lattice. One challenge we foresee is choosing and 
correctly normalizing the operators to use on the lattice. Another challenge is dealing 
with (incorrect or divergent) short- distance behavior of the correlators. We will not 
discuss the issue of overcoming fluctuations to achieve good statistics; instead we 
hope that existing techniques [26] will prove sufficient. 

In general, an operator written in terms of lattice variables will not correspond to 
the continuum operator of interest, but will renormalize and mix with all operators 
with the same symmetry properties. For instance, a proposed lattice implementation 
of T xy will generically be expressed in terms of the true T xy as 

TZ = Z T TZ tin + J2c n O:y, (4.1) 

n 

where O xy are all other operators with the same symmetries as T xy under the lattice 
symmetry group, and Z T , c n are some coefficients. Generally the operators O n 
are higher dimension than T xy and so the c n will carry positive powers of the lattice 
spacing. Therefore, to the extent that we can take the continuum limit the O n should 
be harmless except that they can introduce short-range contributions to correlators. 
However, both the operation of vacuum subtraction and the small momentum limits 
associated with any lattice implementation of &^ G{k)\j Sz ^Q tend to remove sensitivity 
to short distance contributions to the correlators, so we expect this issue to be under 
control. 6 The problem is the renormalization constant Z T , which in general must be 
determined nonperturbatively. 

To evaluate Z T it is useful to recall the physical interpretation of the stress tensor. 
If we make a small change to the geometry, changing g^ u = rj^ to g^ u = rj^ + h^, 
the action changes from 

S g =, = J d A x £ to S g=v+h = J d 4 x (c - + 0(h 2 )^j , (4.2) 

where Cq is the Lagrangian density evaluated assuming = 0. For instance, if we 
modify the lattice action such that the lattice spacing in the x-direction increases, the 
change in the action, to leading order, is — \h xx T xx . Similarly, if —h xy T xy is added to 
the action, the geometry becomes skewed such that the separation between the point 
(0, 0, 0, 0) and the point (0, x, y, 0) is no longer ^ x 2 + y 2 but is \/ x 2 + 2xyh xy + y 2 . 

6 The short distance behavior of the stress tensor two-point function is (T xy {x)T xy (0)) ~ a; -8 . 
For a dimcnsion-6 operator O n , the correlator is (T xv (x)O n (0)) ~ a 2 x~ 1Q . The vacuum subtracted 
value at short distances is 0(T 4 ) by OPE arguments, see Rcf. [27]; hence (T x y(x)O n (0)) T ~ 
a 2 T 4 x- 6 . The short distance contribution to d%{TO n )(k) is - J x x 2 (T(x)O„ {0)) - J x a 2 T 4 x~ 4 
which is 0(a 2 ) and at worst log UV divergent. Higher dimension contaminants carry more powers 
of {a fx) and also contribute at order a 2 . 
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The general strategy for determining the normalization constant Z T on the stress 
tensor is then to include the proposed stress tensor, with small coefficient (— c/2), 
in the action, and to see how much it changes the effective lattice spacing. For 
instance, for a diagonal component such as T xx , one can measure correlation lengths 
along the x-axis and along other lattice axes, or measure the string tension in the xy 
and yz planes. The change in distance determines h xx , and the relation between the 
proposed T xx and the true one is — (c/2)Tp r o poscd = (—h xx /2)T x r x c (unless the change 
also modifies other axis lengths, in which case the proposed T xx is a mixture of T xx , 
T vy etc). 

This technique has been well developed for diagonal components of the stress 
tensor, see for instance [28-31]. To our knowledge it is not as well developed for the 
off-diagonal components. Unfortunately, all of the Kubo relations we have found, 
specifically Eq. (2.10), Eq. (2.11) and Eq. (2.12), involve correlators of off-diagonal 
components of T^ v . But this is easily fixed by performing rotations in our choice 
of axes. For the case of k, we make a 9 = j rotation in the (x,y)— plane, which 
transforms Eq. (2.10) to 



« = ~ l™ ^ (G x t xx (k) - 2G^\k) + G^(k)) . (4.3) 

Of course Q xx ' xx = qv^vv a ^ vanishing k x ^ y by lattice symmetries, so only one needs 
to be evaluated. 

The correlation function found in [22] involved all off-diagonal stress tensors. 
Arnold et al found an expression involving only T yt (z), T xt (z) [32], and we extend it 
to the nonconformal case in the Appendix, see Eq. (A. 35), which we reproduce here: 

X 3 = 2k*-2 lim 6L.a.Gf' te ' ts (p,g). (4.4) 

p z ,q z — >0 

Our expression for A3 still involves the non-diagonal stress tensor T xt . To re- 
express it in terms of diagonal terms we must perform a 45° rotation between the 
x and time axes. This requires a change in the implementation of the periodic 
boundary conditions in Euclidean space, as illustrated in Figure 4. As the figure 
shows, we typically consider a lattice with principal domain running from t — to 
t = (3; equivalently we can say that we consider the field theory over the whole x, t 
plane, but with an identification map which equates every point (x, t) with the point 
(x, t + /3). Introducing rotated coordinates x' — (x + t)/ \[2 and t' — (t — x)j \/2, the 
identification map equates a point at (x',t') with a point at (x f +/3/y/2, tf+/3/ v2). 
We then choose to work on a lattice grid with principal axes along the x' and t' 
directions. We have to pick a principal domain, that is, a region of the plane holding 
exactly one copy of each equivalence class of points under the periodic identification. 
One choice is to stick with the band of points with t G [0,/3); as illustrated in the 
figure, the points labeled 1 are identified, as are the points labeled 2, 3 etc. This 
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Starting Lattice Configuration 
2 3 4 5 6 7 



performing a > 
6=71/4 
rotation 



1 2 3 4 5 6 7 





Figure 3. How to handle a rotation which mixes a time and a space direction, on the 
lattice. 

choice is to consider these points as the boundary points on the lattice, which are 
identified with each other. But the points labeled 2' are also identified, as are the 
points labeled 3' etc. Another choice, as illustrated in the figure, is to choose the band 
with t' G [0,(3/y/2), which has these points as the periodically identified boundary. 
The identification map relates a point on the bottom edge of this band with a point 
on the top edge, but shifted over a distance (3 / a/2 in the x' direction. Note that this 
band is also narrower than the original band; the inverse temperature /3 corresponds 
to the space separation of the identified points, not the extent of the new "time" 
coordinate t! . 

Therefore, implementing lattice gauge theory on a space where the periodic iden- 
tification has a spatial shift corresponds to choosing axes which lie at an angle with 
respect to the (x, t) axes. In this way we can perform the required (x, t) rotation to 
make the stress tensor operators needed in the evaluation of A3 correspond to diag- 
onal components of the stress tensor. Specifically, in terms of the x', t' coordinates 
after this final rotation, A3 is determined by 

A 3 = 2k* - \ p fim q J^- {G tt E ' ttM (p, q) ~ 2G%^(p, q) + G^(p, q)) . (4.5) 
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In a completely analogous way, we find that 

* = « - 2K> - 7 So ^1 ( G| "'" (p - ?) - G|A " (P - ?) ) ' (4 ' 6) 

5 Discussion 

Our central results are presented in Eq. (3.38) and Eq. (3.39). The thermodynamic 
coefficients, unlike the entropy-generating coefficients n etc, do not diverge in the 
weak coupling limit, but remain finite. Both k and A3 are in general nonzero. In 
particular, the previous observation that A3 vanishes in strongly-coupled A^=4 Super- 
Yang-Mills theory [14, 35] appears to be an accident. 

Curiously, our results for k and A 3 , Eq. (3.38) and Eq. (3.39), yield zero when we 
insert the matter content of A/"=4 Super- Yang- Mills theory: iVo = 3N\ and Ni = 4N\. 
Therefore, both coefficients vanish in the weak-coupling limit. As we just mentioned, 
A3 also vanishes in this theory in the strong coupling limit. The fact that A3 vanishes 
in both limits is suggestive that it is strictly zero, but this is not the case; it has been 
shown [36] that A3 is nonzero at subleading order in the large-coupling expansion. 
And it is easy to find other examples of conformal theories where A3 is nonzero. 
For instance, SU(iV c ) gauge theories with fundamental vectorlike matter and with 
the number of flavors Nf slightly below ^N c have weakly coupled, conformal fixed 
points [37]. Since such theories are weakly coupled, Eq. (3.39) applies. Both terms 
are of the same sign, so A3 is certainly not zero for these conformal gauge theories. 
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A Non-conformal hydrodynamics 

In this section we derive Kubo formulae for all the second-order, thermodynamical 
transport coefficients of all non-conformal hydrodynamics. We work in Minkowski 
space and analytically continue to Euclidean space. Consider a hydrodynamic system 
in equilibrium with some arbitrary background of form 

ds 2 = b iliV dx^dx v + hoo(x)dt 2 + hoi(x)dtdx l + hij(x)dx' l dx^ . (A.l) 

In this curved background the expectation value of the stress tensor can be expanded 
about flat space in powers of h^: 

(Tg') h = G% + j d i xG a / flu (x,0)^^- + j d A xd A yGf^ u {x, y, 0) 

(A.2) 
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where the Euclidean Green functions are defined in Eq. (2.8). The corresponding 
Minkowski space expression is 

(Tr) h = G^- Jd*xG:^(x,0)^^ + J d*xd*yG^(x, y, 0) M»)My) ; 

(A.3) 

where here the retarded Green functions are the correlation functions of a T r with 
one or two T a , as is explained in detail in [22]. We can shift from each of these two 
signatures to the other accordingly by noticing that t indices of the Green function 
are multiplied by a factor of % and also for each a index we get an extra minus sign. 

The following expressions for the curvature tensors will be quite handy for the 
derivations of Kubo formulae. To first and second order of perturbations in the 
metric of the form, g^ v (x) = rj^ + h^ v (x), we have (borrowing the result from [38]) 

R a p~/8 = Sfrfapy] — h\^Sx5],p + 'Sa[j,fiS ' p 5],j3 + C(^L) (A-4) 



2 

-^f v {hns,pv — hps^v — h pUt/3 s + h/3 Vtfl s) — (h^ — -I 

.3 



-hrhwfr + 0{h%,) 



3 



and in the above expressions, we have S\s,p = h\s t p + hp\ t s — hps,\. 

We will write thermodynamic variables w M , e, P in an expansion about h = 0; 
^(x) = u^(x) + u^ix) + u^ 2 {x) and similarly for e, P (the barred variables are the 
h = values). It is important to understand the role of the rest frame; w M need 
not equal = (1,0,0,0); rather we must determine by solving conservation 
equations, V = 0, consistently and truncating the expansion. Finally the fluid 
vector in equilibrium must satisfy both a^ v = and V ■ u — 0. 

We also find it useful to take the trace of the energy-momentum tensor in the 
case of non-conformal transport coefficients, 

TH = p(s-\)+3U. (A.5) 



To find the pressure in terms of the background source, we use equations of motion, 
in equilibrium we have 

= V V V^ V (A.6) 

= v,v^ al + v.vx* + v,v M (A^n) 
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where for the ideal fluid, we have 

V,V^ al = uWVvV^e + P) + (e + P)V„V>*V) (A.7) 
+V At (e + P)V v {vru v ) + V„(e + P)V„(u'V) + DP . 

Since it^, — )■ for — )• 0, the first term is identically zero, for the second term 
we get 

(e+PjV^V^iiY) = (e + P) (V„u v VX + «7,V P «" + V,ii"V/ + u"V„VX) 
~ (e + P) {(K x u x ) 2 + r£ Q r>V + + 0(w, ^) 

and we used [V M , V„]w p = R p u". Similarly, for the third and fourth term in 
Eq. (A.7), we get 



V M (e + P)V u (uV) + V„(e + P) V^A") ~ 2<9 a (e + P)r« + ^) . (A.* 



Adding up previous results, finally the pressure reads 

p = P _ + ^gg!fe + p + n) + x _ n + 0(M| ft , , (A . 9) 



where □ = Y^=i&i an< ^ 

x = (e + p) ((r> A ) 2 + r^r^«V + p,X^) + 2 <9 a (e + P)r^V . (a.io) 

This is the generalization of the result in [32]. As pointed out in the main text, some 
Kubo formulae don't directly relate a transport coefficient to the zero frequency and 
momentum limit of Green's functions but they mix these parameters. Throughout 
the next section we'll try to find the simplest setup that can give rise to Green's 
functions manageable for lattice calculations. 

A.l Kubo relation for k and £ 5 

We will start the calculation with £5, by turning on an h xy (x,y) perturbation and 
evaluating (T tt ). But this is e, the energy density of the fluid, we can find it at 
different orders of h^ v from Eq. (A. 9) or by solving the equations of motion directly. 
For illustrative reasons, we do the second approach. Since only h xy (x,y) is nonzero 
we can assume is also a function only of x, y. The viscous tensor is 

T xx _ ld 2 h xy (x,y) c& 2rj ( du\ du x h \ ( ' du\ du x h 



ld 2 h xy (x,y) 


3 


dxdy 


1 d 2 h xy (x,y) 


3 


dxdy 


1 d 2 h xy (x,y) 


3 


dxdy 


-v 1 


(du y h ! < 
^ dx 


-ci 


( du l , « 



n- = ; ^k^(«+«6) + f ^-Vj-Ha«f + * |,(A - n) 

3 \ dx dy J \ dy dx 

n» = -" tt' *' (-2« + ee 5 ) + ^ + ffi )- c ( d -f + dul 



3 \ dx dy J \ dy dx 

f M") , n M = -77^ , IF 2 = -77^ 

ch/ / dx dy 

ax / dxdy 
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here and through the following sections we neglect dissipative terms since they will 
be proportional to the time derivatives of h^ v or fluid vector in general. 

Solving the equations of motion for V M T^ = 0, V^T™ = 0, V M T^ = 0, 
V = 0, we get accordingly, 



dx ^ dt ^ ^ dydx 2 1 3 ^ ^dxdy dx 2 



dP h du x h - d 3 h xy (k \ / d 2 u v h d 2 u x 



3 dx 2 ' dy 2 ' 3 dxdy ) 



Ad 2 u x h d 2 u x h ld 2 u v h 



<9y dt ^ ' dxdy 2 V3 / \dxdy dy 

\3 <9y 2 9a; 2 3 dxdy J 

£<"">-".(£ + 3)-° 

Since we are interested in the zero frequency limit all time derivatives are zero. As 
we can see in the above equations terms proportional to the metric perturbation 
appear, which act as a source for P and for vF~. Higher order terms that include the 
interaction of the fluid vector with background perturbation have been neglected. 
From the first two equations we obtain 

and we know that pressure and energy density are related through P = c 2 e. Similarly 
we get 

(A.14) 

Finally, from Eq. (A. 3), we have 

^ (f + 2^ 5 ) |^ = -\ J *xG°S#(x, 0)M*) • (A-15) 

Fourier transforming and taking the variation of both sides with respect to h xy , we 
get 

^ = "I k l r n ^fmr G ^ "(*0 - i (A.16) 

2 fex.fey-^o ok x dky o 



and we know that 



lim G^' tt (fc , k) = +Gf '"(A*, = 0, k) (A.17) 

feo-)-0 
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So finally in Euclidean space we have 

Similarly the Kubo relation for k in terms of off-diagonal components of stress- 
tensors is given by k = —d\ G x ^: xy {k). In terms of the Euclidean correlator this is 

A. 2 Kubo relation for n* and £ 6 

To find a Kubo relation for £ 6 , we use the perturbation h tt (z). This choice, shifts the 
local rest frame by u l = 1 + l/2h u + 0(h 2 ). Expanding the fluid vector in terms of 
metric perturbation h^, we find the following viscous tensors, 

IF* = IF = (-2k* + k - 3£ 6 + 6^) 

„ d 2 h tt ^ d 2 h tt 

If we assume that the hydrodynamic waves are only functions of z then for V M T MZ = 
we have 

dP h dui , -. Idhu, — i d 3 h tt ( f 6 2k* k \ 

= iF + ir (f+p) -2^ i(f+p) + ^(-! + --3 +& ) (A20 > 

and once again the last term acts as a source for pressure. For (T u ) we have 

The first term in the above relation is a pure gauge. For the second term from 
linear-response we have 

The corresponding Kubo formula for £g will be 

and accordingly in Euclidean space using lim^^o G t ^ tt (k t , k) = —G t ^ tt (k t = 0,k), it 
can be rewritten as 

e, = 26+ £_|_fta *<*•<*). (A.24) 
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This formula determines a linear combination of ^ an d k*. To get both coefficients 
separately we need to look for another relation for k*. We do so by investigating an 
off-diagonal component of energy-momentum tensor. We have 

(T xy ) = (e + P)u x u v + Pg xy + U xy . (A.25) 

Since k* is a coefficient involving the curvature tensor, which first arises at linear 
order in h, we need only work to this order, in which case the first term is zero; and 
if we choose h xy nonzero then Pg xy and A xy U are also zero. Therefore we consider 
hu{x,y)- Then the only contribution comes from w xy , after expansion in the orders 
of h u and u^{x) = u^(x) + u^(x) + 0(h 2 ), we find 

= (V _ *) (A.26) 



dxdy V 2 

If we use Eq. (A. 3) and Fourier transforming, we can write the Kubo formula in 
Euclidean space as 

*" = h\^3&K Gr{k) - (A ' 27) 

Eq. (A. 27) and Eq. (A. 18) involve the same Green function, so we find a relation 
between £ 5 and k*, specifically 

.2,.* , K ( „2 1 



£5 = ~cy + ^ [4 - 3J ■ (A.28) 
A. 3 Kubo relation for A 3 and £ 3 

Now we turn to nonlinear transport coefficients, where we must work to second order 
in h. We begin with A3, which is the traceless contribution arising at second order 
in vorticity. Vorticity is generated by a nonvanishing value of huj] specifically the 
vorticity term for which A3 is a parameter (see Eq. (2.7)) is 

n {i x Q^ x = 1 (6 ij 6 mn - 35 l rn 5l) e mkl e nrs d k h lt d r h st . (A.29) 

The easiest way to proceed [22, 32] is to consider an off-diagonal component of T, 
such as T xy \ then complications involving the pressure (such as those of the last 
subsections) do not arise. However as we have discussed it is most convenient on 
the lattice to use a relation involving a diagonal component of the stress tensor. 
Therefore we will instead consider the vorticity- related contributions to T xx : 

(T xx ) = (e + p)u x u x + pg xx + U xx . (A.30) 

Since we have to keep all terms to the order of 0(h 2 ), we need to know u%, P^, 
and ITff. We will consider nonvanishing h ty (z), which is general enough for U xx to 
contain a A3 dependent term. For this choice u x h vanishes. To find the contribution 
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of P h 2, we look into Eq. (A. 9). Since T^ t is zero for our metric perturbation, the only 
possible contributions come from d il d v ix yw and II. All derivatives other than d z are 
zero, so finally we have 

p h2 = — 7r zz - n = -Tl zz . (A.31) 
Calculating II 22 , U xx and recalling that = 0, the final result reads 

n 22 - h 92hty (V + K } + ( dh °yY( 3 6> A 3 £ 6 K * ,&\, n(h ^ 



Now for Eq. (A. 30) we can write (T xx ) = P + Ii xx = P — II 22 + U xx , which reads 



- P = + ( S -t) 2 (4 + T 1 + 0(ft3) ' (A - 32) 



For the response of the three-point Green's function, we have 

32a /ql \ 2 



(A.33) 

which results in the following Kubo formula for A 3 : 

A 3 = 2*;* + 2 lim _^_^' M (p, g ) , (A.34) 
p*,?«-)-o op z oq z 

and in Euclidean space through \im poqa ^ G t ^ y,xx (p,q) = —G t ^' ty ' xx (p,q), we get 

A 3 = 2*;*- 2 lim ^- G f ty > xx {p,q). (A.35) 
Pz,qz^o op z oq z 

In the main text we use this expression but with x <H> y. 

It's quite easy to find a Kubo relation for £ 3 now, since we have a relation for 
A3. From the trace of stress-tensor we have 

T» = P ( 3 - \ J + 3IT . (A.36) 



Using Eq. (A.31) and the expressions for II 22 and II that can be calculated as 

n = - 2 ^i^r6 + \ (^t) 2 fa ~ 3 ^ + 6) + <^ 3 ) - ( A -37) 

and from the response of the retarded Green's function and then analytically con- 
tinuing to Euclidean space, we finally obtain the Kubo relation for £ 3 : 

6 = 3& - «, + (4 - I) - + c2 f lim o g^-G'n, ,) ■ (A.38) 
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A. 4 Kubo relation for A 4 and £ 4 

The coefficients A4 and £4 arise when the entropy, and therefore temperature, vary 
in space. For this to occur in equilibrium, the gravitational potential h tt must vary 
in space. Therefore we consider perturbations of h tt {x,y) to second order. The 
off-diagonal component of the stress-tensor, (T xy ), is 

xy xy k d 2 h tt dhtt , ( A 4 k K*\dh tt dh tt 

* v = V - ^hu—-^ + k h tt —— : + - + -5- -5- (A.39) 



2 9x9?/ 9x9?/ \^ c s ^ 2 J dx dy 

that after analytic continuation reduces to 

A 4 = -| hm ^-G tt / t ' XV ip,q)-2K* + K. (A.40) 
2 px,?i/^-o op x oq y 

To find £4 we evaluate (T^) for a perturbation of h u {z). ^From Eq. (A. 5), we 
find the pressure through Eq. (A.9). We have T\ z = — l/2d z h tt , T z tt = —l/2d z h, and 
Rtt = -l/2d 2 z h tt - 1/A(d z h tt ) 2 , which result in 

from Eq. (A. 21), and 



X = Xh 



e + P fdhuV d(e + P) h dh t 



dz I dz dz 



_ & + 3^+2k*-k + X 4 \ f dh tt 



4cf 4 6 6c!/ V dz 

£ 6 2^* - k \ d 2 h t 



2 3 J 11 dz 2 ' 
After Fourier transforming and straightforward simplifications, the pressure reads 
p e + P ,e + Pf 1\ p z q z (2K* 2A 4 £ 6 ^ ^ £ 4 ' 



"8 



Inserting the above expressions for P and II from Eq. (A. 42), finally we get 

U = -f lim n -^G^fo g) + 6c 6 s ( | - £ 5 ) (A.44) 
2 p z ,^->o dp z dq z \ 2 J 

L\ ( I ^1 _ & _ ^ * & 

c2j\ + 4V~ + 3rf 2 3 + ^ 5+ rf 
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